A Case Study in Porting Strategy: The Axisymmetric Damage Model

The axisymmetric damage model (ADM) is one of four micromechanical stress field models developed and maintained by N. J. Pagano's group at the AFRL Materials Directorate, Wright-Patterson AFB. These models serve to establish rigorous field theories to compute accurate stresses and energy release rates for fracture mechanics problems in composite bodies. In turn, these solutions, coupled with experimental research, can establish appropriate failure criteria to be used in the design of composite structures for air and space vehicles. Each of these models is based on Reissner's variational theorem (Reissner, 1950). ADM (Pagano, 1991) is associated with axisymmetric cracking and boundary conditions. Others support fundamentally different cases, i.e. axisymmetric problems incorporating the effects of friction (Tandon and Pagano, 1996), transverse and shear loading and associated cracking patterns in composite cylinders (Pagano and Tandon, 1994), and transverse cracking and delamination problems (Schoeppner and Pagano, 1998), but all share the same basic data structures and programming logic which implies that modifications made to ADM are in fact applicable to the other codes as well. ADM eigenanalyses demand a floating point system that supports an unusually wide range of exponent values. The combination of very fine discretizations and the method's requirement for full field analyses can otherwise produce unacceptable levels of error. As a result, production runs have been largely restricted to Hewlett-Packard platforms that provide an exponent range of -4931 to +4932, or about twice that of Cray's floating point format. This unusual characteristic spawned a review of eigenvalue sensitivity and matrix conditioning techniques. Discretizations are based on sets of concentric cylinders that extend into the matrix material and surround an individual core. Each shell has constant thickness and constant inner and outer radii. Damage may be introduced in the form of annular cracks within the constituents and/or debonds between them. The number of shells and their respective thicknesses are dictated by the details of the composite domain, the damage itself, the desired accuracy, and choice of boundary conditions. Cylinders are further subdivided into sections of arbitrary constant lengths, each of which may include different types of damage. Sections may be analyzed independently and then integrated into the global solution downstream of the analysis kernel. Even though this dimension is usually rather small, we chose to implement a multithread design based on it first. The approach held the promise of a significant improvement in performance in exchange for what was viewed as a modest porting effort. The tedious chore of implementing finely grained parallelism was reserved for experimentation with a dependency analysis tool some time later in the development cycle. Key obstacles and the solutions we applied in moving ADM from its largely non-portable, serial form to a multithreaded code are the focus of this report. Performance and scaling characteristics will be covered and a summary of relevant, on-going work will be presented in closing. Introduction Where does one start when faced with the task of modernizing a legacy code? We began with a meeting to review the science and discuss usage and maintenance characteristics. A number of papers describing the implementation changed hands, in particular Brown (1991) and the ADM User Manual (Tandon and Shoeppner, 1992). We made a quick inventory of the source file set (Figure 1). 4 makefiles/executables 76 Fortran source files 16384 lines of code 52 occurrences of COMMON/BRN/ZZ,ZZZ 75 GOTOs 2179 numerical labels Figure 1 Inventory of ADM source file set. Approximately 10,000 CPU hours are spent annually running ADM on the S-class Exemplar at SPAWAR. Complete analyses typically require 10-12 hours of CPU time. The code has always been run on Hewlett-Packard platforms as their compilers support an extended precision (128-bit) floating point system with a particularly wide range of exponent values. The combination of very fine discretizations and the method's requirement for full field analyses can otherwise produce unacceptable levels of error. We closed the meeting by confirming that a coarse-grained form of parallelism was already present in the code. ADM is used to model multiple forms of damage, each of which is discretized separately. These damage regions are denoted independent sections. Parallelizing this part of the code would require variable scoping and the promotion of arrays where indicated. A directives-based implementation was indicated. This is where we began. Summary of Progress To date, four key improvements have been made to the ADM code. A brief description of each follows. Requests for additional information may be sent to [email protected]. • The makefile arrangement was integrated and experiments with compiler optimizations were performed. We quickly learned that the application of higher optimization levels to the original source codes consistently led to incorrect results. • The static memory allocation scheme used in the main program was modified to support automatic array allocation. The result was a single, dynamic executable limited only by physical memory constraints. Coincidently, the main subroutine was made easier to read and refine. • Single processor performance was enhanced, primarily through source code rewrites, compiler optimizations, and the replacement of a bubble sort routine with a simple heapsort code extracted from Press et al. (1996). An I/O routine charged with checking for the presence of NaNs in intermediate and final results was also replaced with a call to the more efficient isnanf() function available on all IEEE machines. • A coarse-grained form of parallelism was implemented through the use of the directive-driven threads model provided by Exemplar SMP platforms (Parallel Model). Reported timings were generated by the Convex performance analyzer (cxpa) and a lightweight threads timer library developed specifically to eliminate overheads associated with recurrent use of cxpa. ADM was originally developed in 1989, but for the most part, the programming style was reminiscent of an earlier era. By this, we mean that there were no comments, no whitespace, and no use of indentation. Virtually all variables were defined implicitly. Nearly every line was strictly uppercase text. For a taste of this style, refer to Figure 2 which illustrates the original form of the independent section loop. DO 201 II=1,NISEC CALL SCTN(AA(NA1),AA(NA2),AA(NA3),AA(NA9),AA(NA19),AA(NA24), .AA(NA25),AA(NA32),AA(NA34),IA(MA2),IA(MA4),IA(MA5),IA(MA7), .IA(MA14),IA(MA15),IA(MA16),IA(MA17),IA(MA18),IA(MA19),II) CALL SCTN2(AA(NA19),AA(NA20),AA(NA21),AA(NA22),AA(NA24), .AA(NA28),AA(NA33),AA(NA34),AA(NA35),IA(MA5),IA(MA7),IA(MA8), .IA(MA14),IA(MA15),IA(MA16),IA(MA22),K2,II) CALL EIGEN(AA(NA1),AA(NA2),AA(NA3),AA(NA9),AA(NA10),AA(NA11), .AA(NA12),AA(NA13),AA(NA14),AA(NA15),AA(NA16),AA(NA17),AA(NA18), .AA(NA19),AA(NA20),AA(NA21),AA(NA22),AA(NA23),AA(NA24),AA(NA25), .AA(NA26),AA(NA27),AA(NA28),AA(NA29),AA(NA30),AA(NA31),AA(NA32), .AA(NA33),AA(NA34),AA(NA35),IA(MA2),IA(MA4),IA(MA5),IA(MA7), .IA(MA8),IA(MA9),IA(MA10),IA(MA11),IA(MA12),IA(MA13),IA(MA14), .IA(MA15),IA(MA16),IA(MA17),IA(MA18),IA(MA19),IA(MA20),IA(MA21), .IA(MA22),IA(MA23),IA(MA24),IA(MA25),IA(MA26),IA(MA27),IA(MA28), .IA(MA29),IA(MA30),IA(MA31),IA(MA32),IA(MA33),IA(MA34),IA(MA35), .II,IOD) CALL CHKNAN(AA(NA11),AA(NA13),AA(NA14),AA(NA19),IA(MA10),IA(MA11), .IA(MA14),II,1,KK,KJ) 201 CONTINUE Figure 2 Original form of the independent section loop. After rewriting the main program and a small number of key subroutines, we found that compilations using higher levels of optimization (primarily +O2 and +O3) generated executables that produced the same results as the original code.